Exploring chronic and transient tumor hypoxia for predicting the efficacy of hypoxia-activated pro-drugs

Hypoxia, a low level of oxygen in the tissue, arises due to an imbalance between the vascular oxygen supply and oxygen demand by the surrounding cells. Typically, hypoxia is viewed as a negative marker of patients’ survival, because of its implication in the development of aggressive tumors and tumor resistance. Several drugs that specifically target the hypoxic cells have been developed, providing an opportunity for exploiting hypoxia to improve cancer treatment. Here, we consider combinations of hypoxia-activated pro-drugs (HAPs) and two compounds that transiently increase intratumoral hypoxia: a vasodilator and a metabolic sensitizer. To effectively design treatment protocols with multiple compounds we used mathematical micro-pharmacology modeling and determined treatment schedules that take advantage of heterogeneous and dynamically changing oxygenation in tumor tissue. Our model was based on data from murine pancreatic cancers treated with evofosfamide (as a HAP) and either hydralazine (as a vasodilator), or pyruvate (as a metabolic sensitizer). Subsequently, this model was used to identify optimal schedules for different treatment combinations. Our simulations showed that schedules of HAPs with the vasodilator had a bimodal distribution, while HAPs with the sensitizer showed an elongated plateau. All schedules were more successful than HAP monotherapy. The three-compound combination had three local optima, depending on the HAPs clearance from the tissue interstitium, each two-fold more effective than baseline HAP treatment. Our study indicates that the three-compound therapy administered in the defined order will improve cancer response and that designing complex schedules could benefit from the use of mathematical modeling.


Supplementary Methods. Computational implementation of the model
The model is defined on a two-dimensional (2D) domain (Ω = [0,200] × [0,100]µ ! ) divided into 100 × 50 square grids of width ℎ = 2.365 individual cells were segmented from the histology image, as shown in Figure 6c in the main text.All cells are of different shapes and sizes, with areas between 4 and 844 ! .
In order to calculate the interstitial fluid flow, we impose zero velocities on the domain boundaries: top  "#$ = [0,0] and bottom  %#""#& = [0,0], and on cell boundaries:  '()) = [0,0], and the influx flow along the left domain boundary:  *+,-.= [1,0].If the forces  / are known, the force-induced fluid velocity according to the regularized Stokeslets method, is given by Eq.( 6) in the main text, and can be rewritten in the following way: where  is the fluid viscosity.For the known velocities O 4 ,  4 R, we need to calculate the unknown forces   at points   = ( / ,  / ), where [ / ] = [ '()) ,  *+,)-.,  "#$ ,  %#""#& ], that yield that fluid flow.Thus, we solve the following matrix equation for   using an iterative GMRES algorithm (gmres.m in Matlab): The resulting fluid velocity flow is shown in Figure 6d in the main text.The zero-velocity values imposed on the cell boundaries guarantee that the compounds (oxygen, sensitizer, and active drug) carried via this advective transport will not cross the cell boundaries accidentally.
To ensure that the cell interior is not penetrated by the compounds unless they are directly absorbed via cell pseudo-receptors, the numerical implementation of the diffusive transport need to be modified to include only the grid points located in the interstitial space and omit these stencil points that are inside the cells.Therefore, the classical four-point stencil for the diffusion equation: ( a To illustrate this concept, two cases of stencils for nonpenetrable diffusion with some grid points located inside the cells are shown in the figure on the right.The appropriate equations are listed below: To ensure stability of the numerical scheme for diffusion, the following Courant-Friedrich-Lewy (CLF) condition much be satisfied: ∆/ℎ !< 1/4, where  is the diffusion coefficient, ∆ is the time step, and ℎ is the grid width.For our simulations, the time step was chosen to be ∆ = 0.15 × 10 9; minutes to satisfy the CLF condition for all diffusive compounds.
Cellular uptake of the diffusive compounds is calculated from four grid points surrounding each of the points belonging to the given cell.To assure that the compound concentration does not attain a negative value, the cellular uptake rate  must satisfy the following condition: 1 − ∆ > 0. Since the cells are of different sizes, the lethal threshold for active drug is determined by normalizing the absorbed drug level by the cell 's area.
This mathematical model has been implemented in the MATLAB system and the associated code is available at the GitHub depository system (https://github.com/rejniaklab/HAP_schedules).

Non-penetrable diffusion stencil.
Examples of a standard 4-point stencil (right-top corner) and two nonpenetrable stencils (7, 8) for numerical implementation of the diffusion equation.
is thus replaced by the a non-penetrable stencil, where  is the number of stencil points that are located outside of the cells: